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I.  INTRODUCTION 


Various  solid-fuel  ramjet  (SFRJ)  tubular  projectiles  have  been  developed, by  the  U.S. 
Army  Ballistic  Research  Laboratory  (BRL).^*^  One  of  these  projectiles;’a  spin- stabilized 
75mm  version  ^see  Figure  1),  has  been  designed  for  use  as  a  tank  gun  training  round 
(TGTR)  for  the  105mm,  M68  tank  cannon.  'iThe  goal  of  this  development  is  to  demonstrate 
a  TGTR  of  low  dispersion  at  three  kilometers  and  with  a  maximum  (safety)  range  of  eight 
kilometers.  The  concept  of  the  TGTR  is  to  use  the  thrust  (hence  low  drag)  of  the  SFRJ 
projectile  to  obtain  a  ballistic  match  with  low  drag  kinetic  energy  projectiles  up  to  three 
kilometers.  Upon  depletion  of  the  propellant  and  choking  of  the  internal  flow,  the  SFRJ 
will  become  a  high-drag  projectile  with  limited  range.  The  40mm  version  of  the  SFRJ 
projectile  has  eilso  been  tested.^  The  flight  characteristics  of  these  projectiles  were  studied 
in  both  free  flight,  using  radar,  and  flight  thru  an  instrumented  ballistic  range. 

In  order  to  investigate  the  internal  flow  in  an  SFRJ,  several  wind  tunnel  tests  were 
conducted  by  the  BRL,  Launch  and  Flight  Division  using  unfueled,  non-spinning,  full- 
scale  models.^  These  tests  provided  internal  sxirface  pressure  measurements  for  a  variety 
of  configurations  at  a  freestream  Mach  number  of  4.03.  Both  choked  (subsonic  internal 
flow)  and  imchoked  (supersonic  internal  flow)  cases  were  studied.  Subsequently,  Nusca, 
Chakravarthy,  and  Goldberg®  used  a  zonal  computational  fluid  dynamics  (CFD)  code  to 
compute  the  internal  flowfield  (choked  and  unchoked  flow)  for  several  75mm  SFRJ  internal 
configurations.  These  computed  results  were  for  unreacting  flow  and  were  compared  with 
measured  internal  wall  pressure  data.  In  addition,  CFD  has  been  applied  to  the  40mm 
SFRJ  projectile.® 

''Computational  efforts  aimed  at  modeling  the  combustion  in  the  SFRJ  are  essential 
to  future  design  modifications  and  applications  of  the  SFRJ  to  high-velocity,  low  drag, 
anti-tank  projectile  systems^  Vos®-^  has  developed  a  finite  chemical  kinetics,  multistep, 
combustion  model  combined  with  a  solution  algorithm  for  p2irabolic  and  elliptic  partial 
differential  equations  (subsonic  flow).  Models  of  this  type  are  candidates  for  use  in  zonal 
CFD  codes.  This  model  has  been  applied  to  a  combustion  chamber  similar  to  that  of  the 
SFRJ.  Instead  of  solid  fuel,  however,  Vos  considered  the  injection  of  hydrogen  gas  (H2) 
into  an  oxygen  (O2)  stream.  Reported  comparisons  of  this  multistep  model  to  a  simple 
one-step  diffusion  flame  modeF  demonstrated  that  even  the  more  basic  approach  provides 
an  adequate  simulation  of  a  reacting  flow. 

^  Afermagen  and  Yalavianchili  R.J.^  **Fir$t  Diagnostic  Tests  of  a  7Smm  Solid  Fuel  Ramjet  Tubular  Projectile,^ 

ARBRL-MR-OStSS,  U.S.  Army  Ballistic  Research  Laboratory,  Aberdeen  Proving  Ground,  MD,  June  19SS.  (AD  .AI30598) 

^ Mermagen  W.H.,  and  Yalamanchili  R.J,,  **Exyerimental  Tests  of  a  105/75mm  Solid  Fuel  Ramjet  Tubular  Projectile,'* 
ARBRL-MR-3416,  U.S.  Army  Ballistic  Research  Laboratory,  Aberdeen  Proving  Ground,  MD,  December  1984.  (AD  B089766) 

^ Nusca,  M.J.,  and  Oskay,  V.O.,  “Comparison  of  Computational  Analysis  with  Flight  Tests  of  a  40mm  Solid-Fuel  Ramjet 
Projectile,"  U.S.  Army  Ballistic  Research  Laboratory,  Aberdeen  Proving  Ground,  MD,  Technical  Report  in  publication. 

*  Kayser  L.D.,  Yalamanchili  R.J.,  and  Trexler  C.,  “Pressure  Measurements  on  the  Interior  Surface  of  a  75mm  Tubular 
Projectile  at  Mach  4,”  ARBRL-MR-37S5,  U.S.  Army  Ballistic  Research  Laboratory,  Aberdeen  Proving  Ground,  MD,  December 
1988. 

^ Nusca,  M.J.,  Chakravarthy,  S.R.,  and  Goldberg,  U.C.,  “Computational  Capability  for  the  Solid  Fuel  Ramjet  Projectile," 
ARBRL-TR-2958,  U.S.  Army  Ballistic  Research  Laboratory,  Aberdeen  Proving  Ground,  MD,  December,  1988. 

^  Vos,  J.B.,  “Calculating  Turbulent  Reacting  Flows  Using  Finite  Chemical  Kinetics."  AIAA  Journal.  Vol.  25,  No.  10, 
October  1987,  pp.  1985- 1372. 

^  Vos,  J.B.,  “The  Calculation  of  Turbulent  Reacting  Flows  with  a  Combustion  Model  Based  on  Finite  Chemical  Kinetics, " 
PhD.  Dissertation,  Dept,  of  Aerospace  Engineering,  Delft  University,  the  Netherlands,  1987. 
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Gielda*  employed  a  time-accurate  MacCormack  predictor-corrector  scheme  to  solve 
the  tmsteady,  axisymmetric  Navier-Stokes  equations  with  five  additional  species  continuity 
equations.  The  code  was  implemented  on  a  CRAY  X-MP/48  computer  and  required  3  to  6 
CPU  hours  per  solution.  The  assumptions  made  in  this  code  include  supersonic  inflow  and 
outflow  boimdaxy  conditions  (i.e.  no  thermal  choking  and  an  oblique  leading  edge  shock) 
and  a  prescribed  solid  fuel  regression  rate  (.11  cm/sec)  and  fuel  surface  temperature  (2000 
R).  The  prescribed  regression  rate  is  used  to  determine  a  constant  injection  velocity  that 
is  applied  to  the  fueled  wall  surface.  A  ttirbulence  model  is  not  employed.  Results  for  inert 
flow  as  well  as  flow  with  solid  fuel  decomposition  and  gas  phase  combustion  axe  achieved. 
The  application  of  this  model  to  the  ramjet  projectile  is  unclear  for  two  reasons.  First, 
the  combustion  chamber  of  a  ramjet  requires  an  inlet  flow  Mach  number  of  about  .2  to  .3. 
Deceleration  of  the  supersonic  freestream,  in  the  diffuser,  may  be  caused  by  a  leading-edge 
normal  shock.®  The  SFRJ  projectile  is  choked  (i.e.  expelled  leading-edge  normal  shock) 
when  combustion  is  occurring  in  flight,  as  observed  in  shadowgraphs  taken  during  range 
testing.®  Thus  a  supersonic  inflow  boundary  condition  for  the  combustion  chamber,  when 
combustion  is  occurring,  is  not  justifiable.  Secondly,  specification  of  an  injection  velocity 
at  the  fueled  wall  precludes  the  use  of  measured  fuel  regression  rates  for  code  validation. 

The  present  steady-state  combustion  model  assumes  a  heterogeneous  mass-controlled 
reaction  between  solid  fuel  (bound  to  the  lateral  wall  of  the  combustion  chamber)  and 
an  incoming  subsonic  airstream.  Figure  2  summarizes  the  main  features  of  steady-state 
combustion  in  the  SFRJ.  Because  the  solution  is  for  steady-state,  solid  fuel  ignition  and 
depletion  events  are  ignored.  The  reaction  between  fuel  and  oxidizer  (air)  yields  products 
that  are  considered  a  single  specie.  The  reaction  occurs  in  one  step  and  is  irreversible.  Fuel 
and  air  are  fully  consumed  in  the  reaction,  and  they  may  not  coexist  at  any  point  in  the 
flowfield.  In  addition,  the  reaction  is  diffusion  controlled  such  that  the  reaction  rate  is  much 
faster  than  the  rate  of  diffusion.  This  diffusion-flame  combustion  model  is  incorporated 
into  the  overall  flow  solution  by  source  terms  in  the  Navier-Stokes  equations  and  equations 
for  the  conservation  of  chemical  species  mass.  Radiation  effects  due  to  combustion  are 
neglected.  The  model  includes  the  effect  of  axial  rotation  of  the  chamber.  Studies  on  this 
and  similar  models  for  combustion  chzunber  flows  have  been  performed. Diffusion-flame 
models  have  been  demonstrated  to  be  as  valid  for  SFRJ  simulations  as  models  with  finite- 
rate  kinetics.^®  However,  finite-rate  kinetics  combustion  models  require  numerical  values 
of  the  reaction  rates  for  all  species.  These  rates  are  not  readily  available  for  complex  solid 
fuels,  especially  those  with  oxidizer  additives. 

®  Gielda,  T.P.,  “Numerical  Simulation  of  Chemically  Reacting  Solid  Fuel  Ramjet  Internal  Flowfields,  ”  Final  Report  for  U.S. 
Army  Ballistic  Research  Laboratory  Contract  No.  DAAL03-SS~D-000l,  February  19SS. 

^  Hill,  P.  and  Peterson,  C.,  Mechanics  and  Thermodynamics  of  Propulsion,  Addison- bVesley,  1965. 

Khalil,  E.E.,  Spalding,  D.B,,  and  Whitelaui,  J.H.,  “The  Calculation  of  Local  Flow  Properties  in  Two-Dimensional  Fur¬ 
naces,”  Intemaiional_Journal_of_Heat_and_Mass_Transf^,  1975,  Vol.  IS,  pp.  775-791. 

Bradshaw,  P.,  Cebeci,  T.,  and  Whitelaw,  J.H.,  Engineering  Calculation  Methods  for  Turbulent  Flows,  Academic  Press, 
New  York,  1981. 

^^Kuo,  K.K.,  Principles  of  Combustion,  John  Wiley  and  Sons,  Inc.,  New  York,  1986. 

Williams,  F.A.,  Combustion  Theory,  Addison- Wesley  Publishing  Company,  Inc.,  1965. 

Vanka,  S.P.,  “Calculation  of  Axisymmetric,  Turbulent,  Confined  Diffusion  Flames,”  AIAA  Journal.  Vol.  24,  No.  3. 
March  1986,  pp.  462-469. 

Kim,  Y.M.,  and  Chung,  T.J.,  “Finite  Element  Analysis  of  Turbulent  Diffusion  Flames,”  AIAA-88-00S2,  Proceedings  of 
the  26th  AIAA  Aerospace  Sciences  Meeting,  January  11-14,  1988,  Reno,  NV. 

Elands,  P.J.M.,  Korting,  P.A.O.C.,  Dijkstra,  F.,  and  Wijchers,  T.,  “Combustion  of  Polyethylene  in  a  Solid  Fuel  Ramjet  - 
A  Comparison  of  Computational  and  Experimental  Results,”  AIAA-88-SO43,  Proceedings  of  the  26th  AIAA  Aerospace  Sciences 
Meeting,  January  11-14,  1988,  Reno,  NV. 
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For  the  present  study  a  Gauss- Seidel  relaxation  scheme  was  applied  to  the  steady  flow 
Navier-Stokes  equations  with  the  local  flow  velocity  constrained  to  be  smaller  than  the  local 
speed  of  soimd  (i.e.  Mach  number  <  1).  Thus  the  governing  equations  are  elliptic.  The 
stream-function  (0)  and  vorticity  (w)  form  of  the  equations  W8«  used.  In  order  to  clearly 
examine  methods  of  incorporating  the  effects  of  solid  fuel  combustion  in  CFD  codes,  the 
present  subsonic  flow  code  was  chosen  over  the  general  flow  code  described  in  Reference 
5.  The  present  combustion  model  can  be  used  with  more  sophisticated  CFD  codes,  and  a 
wider  range  of  flow  conditions  can  be  investigated.  The  present  code  requires  from  0.5  to 
1.0  CPU  hours  per  solution  (inert  and  reacting  flow,  respectively)  on  a  CRAY  X-MP/48 
computer.  The  code  is  therefore  suitable  for  SFRJ  combustor  design  studies. 

The  governing  equations  are  derived  in  Appendix  A.  The  combustion  model,  turbu¬ 
lence  model,  and  boundary  conditions  are  described  in  the  following  sections.  The  solution 
scheme  is  described  in  Appendix  B. 


11.  COMBUSTION  MODEL 


Figure  2  shows  the  assumed  major  features  of  steady-state  combustion  in  the  SFRJ. 
Most  of  these  features  have  been  measured  or  inferred  from  connected-pipe  tests. Behind 
the  injector  a  region  of  recirculating  flow  is  formed.  This  region  serves  to  bring  oxidizer 
from  the  inlet  to  the  solid  fuel  surface.  The  flame  originates  from  this  region,  and  thus  the 
recirculating  flow  acts  as  a  flameholder.  Tests  have  shown  that  combustion  is  not  sustmned 
without  this  flow.^^  The  flame  may  or  may  not  reach  the  centerline  of  the  combustion 
chamber  at  the  exit  plane.  If  it  does  not,  gas  between  the  centerline  and  the  leading  edge 
of  the  flame  consists  of  air  and  products  that  have  been  accelerated  to  a  velocity  above 
the  inlet  value  (assumed  to  be  smaller  that  the  local  sound  speed).  Between  the  flame  and 
the  combustor  wail,  the  gas  consists  of  fuel  and  products.  The  forweird-facing  step,  at  the 
nozzle  entrance,  retards  the  flow,  but  the  flow  does  not  separate. 

Systems  of  chemical  reactions  are  usually  complex.  A  hydrocarbon  fuel,  for  exeunple, 
may  contain  hxmdreds  of  distinct  chemical  species.  During  combustion,  numerous  short¬ 
lived  intermediate  species  are  formed  in  addition  to  the  final  products  that  include  H2O, 
CO,  and  CO2.  Even  if  it  were  technically  desireable  and  computationally  possible  to  include 
8dl  these  species  eind  reactions  within  the  framework  of  a  calculation,  to  attempt  to  do  so 
would  be  defeated  by  2m  ignorance  of  the  thermodynamic,  transport,  and  chemical-kinetic 
properties  of  the  great  majority  of  the  species.  It  is  possible  to  make  useful  predictions 
of  combustion  processes  by  confining  attention  to  an  idealized  reaction  between  fuel  and 
oxidizer,  neglecting  intermediate  compounds  and  considering  simple  reactions  when  fuel 
and  oxidizer  combine  in  a  fixed  mass  ratio. 

Many  of  the  assumptions  commonly  made  in  combustion  models^^aa  are  also  incorpo¬ 
rated  into  the  present  model.  In  particular,  assumptions  made  in  the  Shvab-Zeldovich 
formulation.'^  The  reacting  fluid  is  treated  as  a  continuum  with  equal  mass  diffusivity 

^'’MaeLaren,  R.O.,  and  HoUman,  A.L.,  “SFRJ  Projectile  Development,"  Final  Report,  CSD-S7J0-FR,  United  Technologies 
Chemical  Systems  Division,  Sunnyvale,  CA.,  S6  June  l9St.  (See  also,  Hussain,  S.M.,  and  HoUman,  A.L.,  “Effects  of  Spin  on 
SFRJ  Projectile  Perfo'nnance,"  Monthly  Progress  Report  Ho.  IS,  Contract  No.  DAAAI5-86-C-0008,  1  Dec  -  31  Dec.,  1986.) 
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(Fj)  properties  for  zJl  species.  The  specific  heat  at  constant  pressure  (cp)  of  each  species  is 
assumed  constant  in  the  gas  phase.  Combustion-generated  turbulence,  the  turbulence  gen¬ 
erated  in  addition  to  that  of  the  isothermal  flow,  is  included  by  using  a  phenomenological, 
one-equation  model  (see  Section  III).  The  reacting  solid  surface  is  considered  energetically 
homogeneous  (i.e.  gas-solid  interface  effects  zure  neglected).  Radiation  energy  transfer  is 
asstimed  to  be  negligible.  A  one-step,  mass-controlled  reaction  is  assumed  for  which  the 
reactant  species  are  fully  consumed.  The  reaction  flame  is  considered  to  be  nonpremixed, 
steady,  nearly  isobaric,  and  of  zero  thickness,  i.e.  a  diffusion  flame.^^ 

The  concept  of  the  diffusion  flame  differs  from  that  of  a  premixed  flame  in  which 
reactants  are  perfectly  mixed  prior  to  chemical  reaction.^^  Since  the  solid  fuel  is  bovmd  to 
the  combustor  wall  and  ablates  into  the  airstream  to  react  with  oxidant,  the  diffusion  flame 
model  more  closely  matches  the  physical  environment  of  the  SFRJ.  However,  the  difficulty 
in  modeling  reacting  turbulent  flows  with  nonpremixed  reactants  lies  in  the  handling  of 
the  mean  chemical  production  terms. 

Consider  an  irreversible,  one-step  chemical  reaction  between  fuel  and  oxidizer.  The 
reaction  maintains  a  fixed  stoichiometric  ratio,  and  the  reaction  rate  is  much  faster  than 
the  rate  of  diffusion  so  that  the  reaction  is  diffusion  controlled.^^  The  products  of  combus¬ 
tion  are  treated  as  a  single  species. 

Fuel  -H  n  Oxidizer  — >  (1  -|-  n)  Products 


Conservation  of  chemical  species  (Equation  A24)  can  be  written  for  each  specie  in  the 
above  re6u;tion.  In  particular,  the  mass  conservation  of  fuel  and  oxidizer  are  expressed  as: 
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Where  i2fu  and  Rox  are  the  reaction  rates  per  unit  volume  euid  m  is  a  mass  fraction. 
For  a  general  multi-step  reaction,  these  reaction  rates  are  difficult  to  accurately  merisure. 
Indeed,  the  uncertainty  of  these  data  for  combustion  models  places  a  severe  restriction 
on  the  accuracy  of  the  entire  computation.  However,  for  the  one-step  reaction  described 
above,  the  reaction  rate  (>f  fuel  (iZfu)  and  oxidizer  (i?ox)  negative  (since  these  reactants 
are  consumed),  and  axe  related  by; 

Rfu  ~  Rox/^ 

Rta  =  —  ilpr/(l  -j-  n) 

Rox/^  —  i?pr/(l  -|-  n) 

The  R  terms  can  be  eliminated  from  Equations  1  and  2  by  dividing  Equation  2  by  n 
and  subtracting  the  resultant  equation  from  Equation  1: 


4 


Assuming  that  fuel  Emd  oxidizer  have  identical  diffusion  properties  (i.e.  that  the  dif¬ 
fusion  coefficients,  Fj  =  /leff/Schmidt  Number,  are  equal),  then  Ffu  =  Fox  =  Ffo.  Equation 
3  can  be  written  in  terms  of  (j>h  =  mfu  —  mox/n  as. 


Assuming  that  F^,  =  Fp,  =  Ffp  and  Fox  =  Fpr  =  Fop,  then  Equation  4  can  be  written  in 
terms  of  </>fp  =  mfu  -t-  »7ipr/(l  +  n)  or  <^op  =  rriox/n  -I-  mpr/(l  -f-  n),  as  well. 

In  terms  of  the  general  form  of  the  conservation  equations  (Equation  A25),  the  equa¬ 
tion  for  (i>  =  TTij  is  replaced  with  three  equations,  where: 


^  ^ox/^»  —  1»  —  (Ffo)effi  —  1)  —  0 

^  -~  <^fp  —  ^fu  "t"  ^pr/(l  4”  ~  ~  (Ffp)eff,  C^  =  1,  dij,  —  0 

<f>  =  <i>op  =  mox/n  +  mpr/(l  -f-  n),  a^,  =  1,  =  (Fop)eff,  =  1,  =  0 


Mixture  fractions  (fh,  ftp,  and  /op)  are  formed  from  linear  combinations  of  the  con¬ 
served  variables  (<f>fo,  <^fp,  and  <^op)  and  are  thus  conserved  variables  as  well. 


</>fo  -  (</>fo)ox  _  <i>h  +  1/n  _  /  n  \  /  1  \ 

(</'fo)fu  -  (</'fo)ox  1  +  1/n  U  -t-  n/  VI  -I-  n/ 


(5) 


/fp  = 


fop  - 


<f>{p  {<t>ip)pT 

{4>ip)tM  -  (<^fp)pr 
*^op  (*^op)pr 
(<^op)ox  “  {4>  op)p 


£ _ <^fp  ~  1/(1  +  n)  /l-fn\  _ /Iv 

)pr  l-l/(l-hn)  “I  n  \n) 


^op-  1/(1  +  n) 
1/n  -  1/(1  -I- n) 


^  =  (n  -I-  (j)op  -  n 


(6) 

(7) 


At  the  stoichiometric  air/fuel  ratio  {rriox/Tnt^i  =  n),  fuel  and  oxidizer  are  consumed  (rrifu  = 
mox  =  0)  and  products  are  formed  (mpr  =  1),  thus  the  mixture  fractions  are: 


(^fo)gt  —  0  T  (*/’fp)8t  ~  1  q-n  ’  ~  1  -I-  n  ’  (®) 

(ffo)gt  =  ^  ’  (/op)st  ~  ^  (9) 

The  flame  location  is  indicated  by  the  flow  contour  along  which  /h  =  l/(l-t-n)  and  mpr  =  1. 
At  this  location  both  the  maximum  flow  temperature  and  a  reduction  in  local  flow  density 
occur.^°  The  mass  fractions  (mfu,  mox,  mpr)  are  decoded  from  the  mixture  fractions  by 
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assuming  that  fuel  and  air  are  completely  consumed  in  the  reaction  and  both  fuel  and  air 
ceumot  coexist  at  any  point  in  the  flowfield.  Thus, 


0  <  /fo  <  (/fo),t  :  mfu  =  0,  mox  =  1  -  (1  +  n)/fo,  mpr  =  1  -  mox 

(/fo)8t  <  /fo  <  1  :  mox  =  0,  mfu  =  {^)  /fo  -  (^)  ,  mpr  =  1  -  mf. 

In  stunmary,  the  present  combustion  model  assumes  a  mass-controlled  reaction  be¬ 
tween  fuel  and  oxidizer,  yielding  reaction  products  that  are  treated  as  a  single  specie.  Fuel 
and  air  may  not  coexist  at  any  point  in  the  flowfield.  The  reaction  occurs  in  one  step  and 
is  irreversible.  The  reaction  rates  are  much  faster  than  the  rate  of  diffusion,  i.e.  diffusion 
controlled  reaction.  The  requirement  of  actual  values  of  the  reaction  rates  has  been  elimi¬ 
nated  by  using  a  single  equation  for  the  conservation  of  combined  fuel  and  oxidizer  mass 
(Equation  3).  Radiation  effects  due  to  combustion  are  ignored.  Solid  fuel  ignition  and 
depletion  events  are  not  treated  since  steady  flow  is  assumed. 


III.  TURBULENCE  MODEL 

Modeling  the  effective  viscosity,  /Xeff  =  /^  +  A^ti  for  a  turbulent,  reacting  flow  has  been 
the  subject  of  much  investigation.^^.ia.ia  Pqj.  some  limited  applications,  the  turbulent 
kinetic  energy  and  dissipation  rate  model,  or  k  —  e  model,  has  been  shown  to  be  applicable 
to  the  reacting  flow  in  a  combustion  chamber.*®  The  turbulent  viscosity  is  defined  as: 

tii  =  =  pta-  P  (10) 

where  is  a  function  of  the  turbulent  Reynolds  number,  Ret  =  pk^^^l/p,  with  I  defined 
as  the  turbulence  length  scale.  The  governing  equations  for  k  and  e  cfin  be  written  in  the 
form  of  Equation  A25,  with: 

<f>  —  k,  =  1,  6^  =  (ri)eff,  1,  d(j)  Sic 

0  —  €,  =  1,  hip  —  (Ft)eff,  1,  dtp  Sc 

where  the  diffusion  coefficients  for  k  and  e  are  defined  as  /ieff/Schmidt  Number.  Equa¬ 
tion  10  implies  an  isotropic  turbulent  viscosity  (homogeneous  flow  without  mean  velocity 
variations).*®  With  diffusion  amd  dissipation  terms  in  the  equations  for  turbulent  kinetic 
energy  and  dissipation  rate,  the  isotropic  assumption  represents  a  limitation  on  the  ac¬ 
curacy  of  flow  calculations,  even  without  combustion.  For  example,  normal  stresses  may 
vary  considerably  throughout  the  flowfield.  These  variations,  and  any  implication  that 
they  may  have  for  the  values  of  mean  velocity,  will  not  be  represented  by  this  turbulence 
model.  In  addition,  the  model  has  been  shown  to  yield  good  predictions  in  a  wide  range 
of  boundary  layer  flows,  but  only  limited  testing  has  been  attempted  in  subsonic  flows 
with  recirculation  and  swirl.  For  flows  with  combustion,  the  difficulty  lies  in  an  adequate 
formulation  for  the  source  terms  Sk  and  5^.*® 

**Sa/foJ,  D.R.,  “Cominition-Generated  Turbulence  in  Practical  Combustors,"  A/AA  Journal  of  Propulsion,  Vol.  4,  iVo.  5, 
Sept-Oct  1988,  pp.  S8S-390. 
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Combustion-generated  turbulence,  the  turbulence  generated  in  addition  to  that  of  the 
isothermzJ  flow,  has  been  widely  investigated  using  experimental  procedvires.  For  enclosed 
air-fuel  mixing/reaction  chambers,  the  available  data  has  been  summarized  by  Bzillal.^® 
Combustion-generated  turbulence  is  strongly  affected  by  flow  density,  pressure  2ind  tem¬ 
perature  gradients,  inlet  blockage  ratio  (step  height  of  the  flameholder  and  the  size  of 
the  subsequent  flow  recirculation  zone),  inlet  flow  momentum,  and  flow  swirl  (combustion 
chamber  rotation).  Experiments  in  these  combustors  have  revealed  that,  relative  to  the 
isothermal  flow,  the  combustion  process  increases  the  strength  ^lnd  decreases  the  size  of 
the  flow  recirculation  zone  behind  the  inlet  step.  As  the  size  of  the  recirculation  zone  de¬ 
creases,  combustion-generated  turbulence  increases  sheirply.  Swirl  is  commonly  employed 
in  combustors  to  improve  flame  stability  and  fuel-air  mixing  and  to  shorten  the  flame 
length.  Experimental  evidence  shows  that  combustion-generated  turbulence  rapidly  in¬ 
creases  with  flow  swirl.  Many  of  these  effects  can  be  included  in  phenomenological  models 
that  correlate  well  with  available  data. 

A  phenomenological  model  for  the  effective  viscosity  in  a  confined,  turbulent,  steady 
diffusion  flame  was  proposed  by  Pun  and  Spalding;^®-^^ 

/ieff  =  (rhu^  +  rh  (^w)^)  .^  ^  (11) 

where  AT  is  a  constant  and  is  chosen  to  insure  that  the  proper  value  of  is  prescribed 
at  the  inlet.  The  parameters  D  and  L  are  the  maximum  diameter  and  length  of  the 
combustion  chamber,  respectively.  The  inlet  axial  and  azimuthal  velocity  are  denoted  u  and 
V,  respectively.  The  inlet  swirl  velocity  is  rv.  This  mode)  assumes  that  ^eff  increases  with 
p,  D,  and  the  momentum  (or  kinetic  energy)  of  the  fluid  entering  the  chamber  (including 
flow  swirl),  but  decreases  with  increasing  chamber  length,  L.  This  model  has  been  studied 
for  simulation  of  steady,  confined,  turbulent  flames.^® 


IV.  BOUNDARY  CONDITIONS 

For  an  SFRJ,  the  present  combustion  model  considers  only  the  flow  within  the  fueled 
combustion  chamber.  The  boundaries  of  this  chamber  are  the  inlet  at  the  injector,  the 
exit  at  the  nozzle  throat,  the  symmetry  axis  about  which  the  chamber  may  rotate,  and 
the  solid  wall  (Figure  2).  Since  the  governing  equations  (Equation  A25)  are  elliptic,  the 
conditions  along  these  boundaries  must  prescribe  values  of  the  dependent  variables,  the 
gradient  of  the  dependent  variables  in  the  normal  direction,  or  an  algebraic  relation  which 
connects  the  values  of  the  dependent  variables  to  the  normal  component  of  velocity.  The 
sensitivity  of  the  flow  solution  to  the  boundary  conditions  has  been  examined. 

At  the  inlet  plane,  radial  profiles  of  all  dependent  variables  —  w,  or  f,h,  and  u, 
-  as  well  as  values  for  T,  p,  p,  and  are  specified.  For  combustion  chcunber  flows,  the 
sensitivity  of  the  flow  solution  to  the  shape  of  the  initial  velocity  and  temperature  profiles 

’®P«n,  W.M.,  and  Spalding,  D.B.,  “A  Proetdart  for  Predicting  the  Velocity  and  Temperature  Distributions  in  a  Confined. 
^967^’  ^*’^**^*"*’  Diffusion  Flame,"  Imperial  College,  Mechanical  Engineering  Department  Report  SF/TN/ll. 

Patankar,  S.V.,  and  Spalding,  D.B.,  “A  Computer  Model  for  Three-Dimensional  Flow  in  Furnaces,"  Proceedings  of  the 
nth  International  Symposium  on  Combustion,  Pennsylvania  State  University,  August  20-25,  1972,  pp. 605-614. 
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has  been  found  to  be  insignificant. It  is  assumed  that  the  fiow  at  the  inlet  plane  consists 
of  air  (/fo  =  0,  /fp  =  —  1/n,  fop  =  1)  and  that  the  combustion  process  does  not  effect  the 
inlet  fiow.  In  addition,  a  subsonic  inlet  fiow  assumption  is  used  in  accordeince  with  the 
elliptic  nature  of  the  governing  equations.  These  assumptions  axe  not  unreEisonable  for 
the  SFRJ  projectile.  The  solid-fuel  ramjet  projectile  is  designed  for  subsonic  combustion^^ 
and  only  limited  evidence  for  the  leakage  of  combustion  products  forward  of  the  injector 
hza  been  observed.*  Hence,  elliptic  solvers  can  be  used. 

The  exit  plane  is  located  at  the  nozzle  throat  where  the  fiow  is  assumed  to  be  subsonic. 
The  streeimlines  at  the  exit  plane  are  assumed  to  be  parallel  to  the  symmetry  axis;  thus 
flow  gradients  along  the  streamlines  are  zero.  These  assumptions  are  reasonable  since 
experience  for  large  Rejmolds  numbers  has  shown  that  the  exact  nature  of  the  exit  plane 
boimdaxy  conditions  has  little  effect  on  the  flowfield  solution  when  convection  is  significant 
(i.e.  large  inlet  mass  flow).^° 


For  mass  continuity,  the  symmetry  axis  is  considered  to  be  a  streamline  of  the  flow, 
thus  xl>  =  constant.  Along  the  symmetry  axis  r  =  0,  thus  d^/dz  =  dtp /dr  =  0  via 
Equations  A1  and  A2.  The  value  of  tp  along  this  boundary  can  be  determined  from  values 
of  tp  adjacent  to  the  boimdary  using  a  one-sided  finite  difference  for  dtp /dr  at  the  zixis. 
From  Equation  A3,  the  boundary  value  for  u  is  zero.  The  aixis  boimdary  values  for  the 
remaining  dependent  variables,  <^,  are  determined  from  d<p/dr  =  0. 


The  no-slip  condition  (u  =  w  =  Oy  and  v/r  =  fl,  where  Q  is  the  wall  spinrate)  is  applied 
to  the  solid  walls.  Therefore  tp  =  constant,  via  Equations  A1  and  A2.  For  convenience 
tp  =  0  is  chosen.  One-sided  finite  differences  for  dtp /dr  and  dtpjdz  zxe  used  in  Equation 
A3  to  determine  the  wall  value  for  w.  For  an  inert  wall  (unfueled),  the  normal  gradient  of 
all  mixture  fractions,  df/dn,  are  set  to  zero.  For  a  fueled  wall,  the  boundeu^^  condition  is 
based  on  the  assumption  that  the  wall  is  a  source  of  fuel  (/fo  =  1,  /fp  =  1,  /op  =  — n). 
The  wall  temperature  is  set  to  the  vaporization  temperature  of  the  fuel,  T^aU  =  ^vap-  The 
rate  of  burning  (regression  rate)  on  the  fueled  wall  varies  as  a  function  of  position  along 
the  wall  and  is  determined  from  the  temperature  gradient  normal  to  the  w^,^® 


^wall 


(12) 


where  is  the  local  heat  flux  to  the  wall  (-«5T/5r)waii,  />waU  is  the  density  of  the  fueled 
wzdl,  2ind  hy  is  the  heat  of  vaporization  of  a  unit  mass  of  solid  fuel.  Values  of  the  thermal 
conductivity  (/c),  density,  and  heat  of  vaporization  can  be  determined  for  most  solid  fuels. 


V.  VALIDATION 


The  present  steady-state,  elliptic  flow  solver  and  combustion  model  can  be  validated 
using  available  data  for  ramjet-type  configurations  (e.g.  sudden-expansion  combustors). 
Detailed  experimental  data  for  reacting  internal  flows  in  these  configurations  are  rare,  and 
the  current  simulation  is  aimed  at  overall  flow  properties.  Comparison  with  measured  data 
focuses  on  three  features:  injector  flows,  fuel  regression  rate,  and  thrust. 

The  reaxward-facing  step  at  the  inlet  (i.e.  the  injector)  is  required  to  initiate  and 
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sustain  the  combustion  flame  in  am  SFRJ.^’^  Within  this  zone,  air  is  mixed  with  fuel,  and 
the  solid-fuel  regression  rate  is  a  maximum.  Therefore,  proper  simulation  of  the  flow 
reattawdiment  point  on  the  fueled  wall  is  essential. 

The  solid  fuel  regression  rate  (r)  in  SFRJ  combustors  has  been  measured  using  several 
methods.  Equation  12  indicates  that  r  is  proportional  to  the  temperature  gradient  normal 
to  the  fueled  wall.  Accturate  simulation  of  this  quantity  is  essential  to  SFRJ  design  studies 
where  the  fuel  bum  time  must  be  selected  to  achieve  required  flight  characteristics. 

Prediction  of  the  thrust  generated  by  a  SFRJ  is  also  essential  to  design  studies  and 
will  be  addressed  in  Section  VI.6. 


1.  Reattachment  Length  Data 

Lee*^  conducted  a  computational  study  of  recirculating  laminar  and  turbulent  flows 
through  an  abrupt  pipe  expansion.  The  flow  was  chemically  inert.  This  study  employed  the 
stream  function- vorticity  equations  in  flnite  difference  form  and  a.  k—e  model  for  turbulence 
closure.  A  Gauss-Seidel  solution  scheme  was  used.  Results  were  generated  for  a  variety  of 
inlet  flow  profiles  including  uniform  or  plug  flow  and  fully  developed  flow.  Flow  Reynolds 
numbers  (based  on  the  inlet  pipe  diameter),  ranged  from  10  to  200  for  laminair  boundary 
layers  and  1000  to  1  million  for  turbulent  flows.  The  ratio  of  internal  pipe  diameter  to  inlet 
diameter  was  2  (step  height  equivalent  to  25%  of  the  combustor  diameter).  Comparison 
with  experimental  data  from  several  soiurces  was  used  to  validate  the  code.  The  present 
code,  without  combustion  model,  was  compared  to  Lee’s  calculations  of  flow  reattachment 
length  for  laminar  Reynolds  numbers  and  uniform  inlet  flow.  Turbulent  calculations  were 
not  compared  since  a.  k  —  e  turbulence  model  was  not  employed  in  the  present  code.  Figure 
3  shows  a  nearly  linear  variation  of  reattachment  length  with  Reynolds  number.  Agreement 
between  two  independent  calculations  is  good. 

So  and  Ahmed^^  conducted  an  experimented  investigation  of  inert  flow  in  a  sudden- 
expansion  combustor.  Low  speed  flow  (13.1  ft/sec)  in  a  long  tube  with  a  well  contoured 
inlet  and  nozzle  was  meeisured  with  a  laser-Doppler  velocimeter  (LDV)  for  various  values 
of  combustor  axial  spinrate.  The  ratio  of  combustor  diameter  to  inlet  diameter  (at  the  end 
of  the  inlet  nozzle)  was  1.5  (step  height  equivalent  to  16%  of  the  combustor  diameter).  At 
the  entrance  to  the  combustor,  the  flow  had  been  accelerated  to  34.8  ft /sec.  The  Reynolds 
number  based  on  this  velocity  and  the  inlet  diameter  was  31,000.  Figure  4  indicates 
that  the  flow  recirculation  zone  is  compressed  by  combustor  spin.  The  flow  reattachment 
length  (measured  from  the  step)  is  normalized  by  the  step  height  above  the  combustor 
wall.  Comparison  of  this  data  with  results  from  the  present  code  axe  also  shown.  While 
the  zero  spin  value  is  reproduced  by  the  present  code,  calculations  show  that  the  variation 
with  spinrate  is  less  signiflcemt  than  indicated  by  the  data. 

Let,  J.J.,  “Heat  and  Mate  Transfer  Calealationt  for  Reeirealating  Laminar  and  Tarialent  Flows  in  an  Abrapt  Pipe 
Expansion,"  AIAA-88-3791,  Proceedings  of  the  1st  AJAA/ASME/SJAM/APS  National  Fluid  Dynamics  Congress,  July  S5-S8, 
1988,  Cincinnati,  OH. 

So,  R.M.C.,  and  Ahmed,  S.A.,  “Rotation  Effects  on  Axisymmetric  Sudden- Expansion  Flows," 

AIAA  Journal  of  Propulsion  and  Power,  Vol.  d.  No.  3,  May-June  1988,  pp.  S70-S7S. 
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2.  Regression  Rate  Data 


Holzman  and  Husian^^  at  United  Technologies  Chemical  Systems  Division  (UTC- 
CSD)  have  measured  the  solid  fuel  regression  rate  in  an  SFRJ  using  a  connected- pipe  test 
stand.  The  internal  geometry  of  the  model  closely  resembled  that  of  the  75mm  SFRJ 
(Figure  1).  The  solid  fuel  used  in  the  tests  was  the  same  as  that  used  in  the  SFRJ 
projectile.  A  typical  irdet  mass  flux  for  the  tests  was  1814  grams/sec  (4  Ib/sec)  at  an  inlet 
velocity  of  approximately  9.1  m/sec  (30  ft/sec)  and  a  temperature  of  1840  R.  The  ratio 
of  the  combustor  diameter  to  the  inlet  diameter  was  1.42  (step  height  equivalent  to  15% 
of  the  combustor  diameter).  Measurements  were  taken  for  combustor  axial  spinrates  from 
0  to  35,000  rpm  (3665  rad/sec).  The  solid  fuel  regression  rate  was  determined  from  the 
meeisured  fuel  grain  regression  and  the  test  run  time. 

The  solid  fuel  consists  of  four  main  components:  30%  ammonium  perchloride  (NH4CIO4), 
28.5%  hydroxl- terminated  polybutadiene  (C4H6.052O0.005),  40%  Escorez  (C7.398Hn.059),  and 
1.5%  Catocene  (C27H32Fe2).  The  main  products  of  combustion  are  hydrocarbon,  oxygen, 
water  and  nitrogen. For  this  study  the  assumed  products  are  67%  CH2,  16.3%  O2, 
12.1%  H2O,  and  4.9%  N2.  Thus,  the  molecular  weights  of  the  fuel,  oxidizer  and  products 
are  97.7,  28.8,  and  18.2  grams/mole,  respectively.  The  fuel  density,  heat  of  vaporization,  va¬ 
porization  temperature,  and  stoichiometric  air/fuel  ratio  are  1.189  grams/cc  (74.2  Ib/ft®), 
2.8998  X  10^  Joules/kg  (12465.6  BTU/lb),  1755  R,  and  9.152  respectively,  as  determined 
by  UTC-CSD. 

Figure  5  shows  the  fuel  regression  rate  distribution  along  the  fueled  wall.  The  mea¬ 
sured  data  are  for  chamber  spinrates  between  0  and  35,000  rpm.  In  this  case,  the  solid  fuel 
did  not  contain  ammonium  perchloride  (AP)  as  an  oxidizer  additive.  The  effect  of  axial 
spin  of  the  combustion  chamber  was  found  to  be  insignificant.  The  present  computations 
for  solid  fuel  without  AP  reproduce  average  or  mean  fuel  regression  rates  reasonably  well. 

At  the  beginning  of  the  fuel  grain  behind  the  injector,  the  computed  regression  rate  is 
larger  than  measured.  However,  the  overall  trend  is  similar  including  the  peak  regression 
rate  at  the  flow  reattachment  point.  There  are  two  possible  explanations  for  this  dis¬ 
crepancy.  In  the  recirculating  flow,  radiation  effects  are  very  significant.^^  However  the 
present  code  ignores  this  effect.  Turbulence  modeling  for  the  recirculating  flow  in  SFRJ 
combustors  is  difficult®  and  may  not  be  properly  modeled  using  the  present  code.  At  the 
end  of  the  combustor,  only  the  trend  of  regression  rate  is  simulated.  In  this  region  the 
fuel  grain  joins  the  nozzle,  producing  a  corner  flow.  Possible  errors  in  the  specified  wall 
boundary  conditions  as  well  as  errors  in  corner  flow  turbulence  modeling  may  explain  the 
discrepancy. 

Elands  et.  al.'®  and  Netzer^**  have  measured  the  solid  fuel  regression  rate  in  an  SFRJ 
chamber  using  a  connected-pipe  test  facility.  Polyethylene  (PE)  fuel  was  used  in  a  small- 
scale  SFRJ  configuration.  The  inlet  diameter  was  15mm  (.6  inches),  and  the  maximum 
chamber  diameter  was  45mm  (1.8  inches)  making  the  ratio  of  combustor  to  inlet  diameter 
3  (step  height  equivalent  to  33%  of  the  combustor  diameter).  Measurements  of  local 

Kubota  N.,  "Survey  of  Rockti  ProptUanii  and  Thetr  CombusUon  Characiertstxcs/'  Fundamentals  of  Solid  Propellant 
Charactertstics.  AIAA  Education  Series,  1984- 

^*Netzer,  A.,  and  Gany,  A.,  "Burning  and  Flameholding  Characteristics  of  a  Afintature  Solid  Fuel  Ramjet  Combustor,” 
AIAA-88-S044f  Proceedings  of  the  26th  AIAA  Aerospace  Sciences  Meeting,  January  U-14,  t988,  Reno,  NV. 
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regression  rate  on  the  fuel  grain  were  taken  using  a  non-intrusive  ultrasonic  method.  The 
effect  of  inlet  air  mass  flow,  inlet  mr  temperature,  and  chamber  pressure  on  fuel  regression 
rate  were  studied. 

For  comparison  to  the  present  code  the  values  of  fuel  molecular  weight,  density,  va¬ 
porization  temperature,  and  heat  of  vaporization  were  taken  as  28.05  grams/mole  (for 
C2H4),  .8008  grams/cc  (50.0  Ib/ft^),  1440  R,  and  4.0009  x  10®  Joules/kg  (1719.9  BTU/lb), 
respectively.  The  inlet  air  temperature  was  540  R.  Figure  6  shows  that  the  average  fuel 
regression  rate,  over  the  fuel  grzun  length,  increases  slightly  for  increased  inlet  air  mass 
flow  rate.  The  results  of  a  finite  rate  chemistry  code,  reported  in  Reference  16,  and  the 
present  diffusion  flame  code  are  also  shown.  The  present  code  does  not  require  thermo¬ 
dynamic  burning  rates  for  the  PE  fuel.  However,  since  PE  is  a  simple  hydrocarbon  fuel, 
these  rates  are  available  for  finite  rate  models.  Both  combustion  models  reproduce  the 
trend  and  magnitude  of  the  data.  The  computational  grid  used  by  the  present  code,  338  x 
44,  was  more  dense  than  the  grid  used  by  the  finite  rate  code,  20  x  20.  The  present  code 
required  0.7  to  1.0  hours  per  solution  on  a  CRAY  X-MP/48  computer. 


VI.  RESULTS  AND  DISCUSSION 


The  present  computational  fluid  dynamics  code  and  combustion  model  was  used  to 
simulate  the  reacting  flowfield  in  the  axisymmetric  combustion  chamber  of  a  75mm  SFRJ 
projectile  (Figure  1).  The  solid  fuel  used  in  this  chamber  is  the  same  as  that  described  in 
Section  V.2.  In  this  section  detailed  computed  results  are  presented  for  the  flame  location, 
flow  velocity  and  mass  fraction  profiles,  and  thrust. 


1.  Inlet  and  Fueled- Wall  Boundary  Conditions 

Several  experimental  and  computational  studies^-'*-®  of  SFRJ  internal  flows,  without 
combustion,  have  revealed  subsonic  flow  in  the  combustion  chamber  for  small  nozzle  throat 
diameters.  Consider  an  inert  SFRJ  with  a  1.1  inch  nozzle  throat  diameter  and  a  freestream 
Mach  number  of  4.03.  The  nozzle  is  choked,  and  a  normal  shock  is  attached  to  the 
leading  edge  of  the  inlet.^-®  These  conditions  are  also  characteristic  of  a  fueled  SFRJ  with 
combustion  but  with  a  larger  nozzle  throat  diameter  (1.6  inches)  and  the  same  Mach 
number.'^  However,  the  cold  flow  measurements  and  computations  for  an  inert  SFRJ  with 
a  1.6  inch  nozzle  throat  diameter  show  supersonic  internal  flow  (unchoked).'*’®  In  a  fueled 
SFRJ  with  combustion,  it  can  be  assumed  that  the  1.6  inch  diameter  nozzle  yields  choked 
internal  flow  (i.e.  leading-edge  normal  shock)  just  as  the  1.1  inch  diameter  nozzle  in  an 
inert  SFRJ.  Hence,  computational  results  for  the  1.1  inch  nozzle  (subsonic  internal  flow 
and  choked  nozzle)®  with  non-reacting  flow  can  be  used  to  provide  inlet  plane  flow  profiles 
for  the  present  combustion  model. 

Computations  are  presented  for  the  75mm  SFRJ  with  1.7  inch  injector  and  1.6  inch 
nozzle  throat  diameters.  Non- reacting  flow  computed  results  (wind  tunnel  test  conditions 
with  Mach  number  4.03  and  Reynolds  Number  20  million  per  foot)  for  a  1.1  inch  nozzle 
diameter  show  that  the  flow  velocity,  pressure,  density  and  temperature  are  nearly  constant 
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across  the  entrance  to  the  combustion  chamber  (except  very  close  to  the  wall).  Average 
viilues  for  these  flow  variables  are  84.40  ni/sec  (276.92  ft/sec),  18192.8  kg/m^  (3725.76 
Ib/ft^),  2.201  kg/m®  (.00427  slug/ft®),  and  510.2  R,  respectively.  These  values  are  taken  as 
sufflcient  upstream  conditions  for  the  present  simulation.  The  upstream  flow  is  assumed 
to  consist  entirely  of  air  and  is  uneffected  by  flow  conditions  downstream.  The  Reynolds 
number  is  3  million  per  foot,  based  on  the  entrance  flow.  A  The  Prandtl  and  Schmidt 
numbers  are  assumed  to  be  unity.  A  Prandtl  number  (Pr)  of  unity  is  considered  adequate 
for  gaseous  flows,  even  with  combustion^^.  A  Schmidt  number  (Sc)  of  unity  implies  that 
the  viscous  transport  of  momentum  acts  in  the  same  manner  as  the  diffusive  transport  of 
mass.  Thus  the  Lewis  number,  Le  =  Sc/Pr,  is  also  unity  which  implies  that  the  energy 
transfer  due  to  conduction  is  equal  to  that  due  to  diffusion.  These  assumptions  are  typical 
of  combustion  models  such  as  the  Shvab-Zeldovich  formulation.^® 

The  mixture  fraction  boundary  conditions  on  the  fueled  wall  are  /fo  =  Ij  /fp  =  1, 
and  fop  =  — n.  Thus,  oxidizer  and  products  are  not  present  at  the  fueled  wall.  The 
wall  temperature  is  equivalent  to  the  fuel  vaporization  temperature,  which  is  constant. 
A  temperature  dependent  fuel  mass  fraction  boundary  condition  on  the  wall  could  be 
implemented  in  the  code.  However,  experimental  evidence  has  shown  that  this  condition 
does  not  represent  the  nature  of  combustion  at  the  fueled  wall^^.  The  combustion  of  fuel 
over  the  entire  wall  is  started  and  sustained  by  the  recirculating  flow  behind  the  injector 
step.  The  entire  length  of  fuel  grain  is  either  burning  or  extinguished  (although  the  burning 
rate  may  be  a  function  of  position  along  the  weill).  Thus,  after  combustion  is  started,  the 
fuel  grain  is  insulated  from  the  temperature  of  the  airstream.  Since  the  present  combustion 
and  flow  model  assume  a  steady  state,  fuel  ignition  and  burnout  (due  to  depletion  of  the 
fuel  supply)  are  not  considered. 


2.  Computational  Grid 

The  computational  grid  used  to  discretize  the  flowfield  consists  of  185  points  between 
the  inlet  plane  (at  the  leading  edge  of  the  injector,  i  =  0)  and  the  nozzle  throat  (x  = 
5.942  inches).  Between  the  centerline  and  the  waU  of  the  combustion  chamber  (1.1  inches 
maximum  radius),  35  points  are  distributed.  The  grid  point  distribution  is  essentially 
uniform  in  both  directions.  Radi2il  wall  spacings  of  .02  inches  at  the  injector,  .03  inches 
at  the  center  of  the  combustion  chamber  (maximum  radius),  and  .02  inches  at  the  nozzle 
throat  are  used.  In  the  streamwise  direction  an  average  grid  spacing  of  .03  inches  is  used. 
The  grid  distribution  should  be  adequate  since  the  turbulence  model  is  not  based  on  wall 
functions,  which  require  exceptionally  fine  wall  grid  clustering.  A  finer  grid  with  essentially 
half  the  wall  spacing  did  not  change  the  results  to  a  noticeable  degree,  but  the  computer 
run  times  where  substantially  increased. 

Converged  solutions  were  obtained  after  about  843  iterations  and  approximately  15 
minutes  of  CPU  time  on  a  CRAY  X-MP/48  computer.  Convergence  was  determined  using 
Equation  B16  with  Ai  set  to  .005.  No  under-relaxation  was  required,  thus  Oyr  =  1- 
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3.  Mean  Flow  Features 


Figure  7  shows  the  velocity  vectors  in  the  vicinity  of  the  injector.  The  flow  directed 
toward  the  fueled  wall  separates  behind  the  bax:kward-facing  step.  The  incoming  flow 
is  deflected  toward  the  symmetry  axis  of  the  chamber.  The  separated  flow  is  decelerated 
near  the  wall  while  the  deflected  airstream  retains  most  of  the  inlet  velocity.  The  separated 
flow  reattaches  on  the  fueled  wall  at  about  .25  inches  from  the  inlet,  thus  the  length  of  the 
recirculation  zone  is  about  .15  inches  or  3%  of  the  chamber  length.  This  is  the  only  region 
of  recirculating  flow  in  the  chamber. 

Figme  8  shows  the  velocity  vectors  in  the  vicinity  of  the  nozzle  entrance.  The  flow 
near  the  fueled  wall  stagnates,  but  it  does  not  separate  in  the  comer  of  the  forward-facing 
step.  The  flow  away  from  the  wall  turns  sharply  toward  the  combustor  centerline.  As  a 
result,  mean  flow  plays  a  significant  role  in  transporting  fuel  from  the  wall  to  the  core  flow. 
In  contrast,  turbulent  flow  convection  is  the  dominant  transport  mechanism  along  most  of 
the  fueled  wall  where  the  flow  is  parallel  to  the  wall.  In  the  recirculation  region  behind 
the  injector,  recirculation  works  against  convection  in  transporting  fuel  into  the  core  flow. 


4.  Combustion  Flame 

Figure  9  shows  the  temperature  contours  across  the  combustion  chamber.  The  highest 
temperature  contoxu:  shown  is  T  =  4297  R  which  is  equivalent  to  the  flame  temperature 
assumed  in  the  computation.  The  concentration  of  high  temperatures  at  the  beginning  of 
the  chamber  indicates  that  the  flame  is  attaxihed  to  the  edge  of  the  injector.  Thus  the  flame 
is  “held”  by  the  injector  or  flameholder.  The  temperature  gradients  at  the  flameholder  are 
large  (i.e.  closely  spaced  contours).  This  result  is  typical  of  diffusion-flame  combustion 
models  where  the  reaction  rate  time  lag  is  zero.  The  flame  extends  to  the  nozzle  and 
intersects  the  chamber  exit  plane  about  midway  between  the  axis  and  the  nozzle  wall.  The 
flow  between  the  flame  and  the  chamber  axis  rapidly  cools.  The  flow  between  the  flame 
and  the  fueled  wall  cools  graduailly  to  the  fuel  vaporization  temperature  of  about  1773  R. 
Most  of  the  heat  produced  by  the  reaction  is  found  between  the  flame  and  the  wall.  These 
temperature  contours  are  functions  of  the  fuel  molecular  weight  and  the  thermodynamic 
data  used  in  the  computation. 


5.  Flow  Profiles 

Figure  10  shows  streamwise  (u)  velocity  profiles  for  four  stations  between  the  injector 
imd  the  nozzle  throat.  The  first  station  (x  =  .136  in.)  is  located  at  about  the  center  of  the 
recirculation  region,  behind  the  injector.  The  second  (x  =  2.038  in.)  is  located  at  about 
chamber  midpoint.  The  third  (x  =  5.572  in.)  is  located  near  the  nozzle  entrance.  The 
fourth  (x  =  5.942  in.)  is  located  at  the  nozzle  throat.  The  no-slip  wall  boundary  condition 
and  the  zero-radial-gradient  axis  boundary  condition  are  illustrated.  Salient  features  of 
the  developing  flow  are:  separation  behind  the  injector  step  (x  =  .136  in.),  flow  stagnation 
(but  no  separation)  near  the  wall  before  the  nozzle  entrance  (x  =  5.572  in.),  and  flow 
acceleration  along  the  chamber  axis  (from  about  310  to  1000  ft/s).  The  radi8d  location 
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of  the  flame  is  indicated  by  a  maximum  in  each  velocity  profile,  which  moves  toward  the 
chamber  axis  as  the  location  of  the  profile  moves  downstream. 

Figures  11  thru  13  show  profiles  of  mass  fractions  (mfo,  mfp,  and  mop)  at  three  axial 
locations.  At  each  eixial  station,  the  flame  location  is  indicated  by  mp,  =  l,mfu  =  0,  and 
mox  =  0.  Since  it  is  assumed  that  the  reaction  goes  to  completion,  the  mass  fraction  of  fuel 
is  zero  between  the  flame  and  the  symmetry  axis,  and  the  mass  fraction  of  oxidizer  is  zero 
between  the  flame  and  the  fueled  wall.  The  flame  location  changes  from  about  y  =  .85  to 
y  =  .50  inches  from  the  axis  between  the  injector  and  the  nozzle  throat. 

Figures  14  and  15  show  the  temperature  and  density  profiles.  At  each  axial  location, 
the  radial  location  of  the  flame  is  indicated  by  the  temperature  peak.  This  peak  is  about 
4300  R  and  is  nearly  constant  from  the  injector  to  the  nozzle  throat.  The  density  profiles 
also  serve  as  a  means  to  track  the  radial  location  of  the  flame.  The  flow  density  reaches  a 
minimum  in  the  flame,  about  7%  of  the  inlet  air  density.  The  flow  density  near  the  fueled 
wall  varies  considerably  along  the  length  of  the  chamber. 


6.  SFRJ  Flight  Performance 

Mermagen  and  Yalamanchili^ have  conducted  free-flight  tests  of  75mm  SFRJ  pro¬ 
jectiles.  A  Hawk  radar/doppler  velocimeter  was  used  to  measure  the  velocity  histories 
of  rounds  with  different  injector,  nozzle  throat  diameters,  and  solid  fuels  (including  inert 
rounds).  All  fueled  rounds  demonstrated  successful  autoignition  with  burn  times  vary¬ 
ing  between  1.6  and  2.5  seconds.  Suitable  repeatabihty  of  fuel  burn  time  and  velocity 
performance  was  reported. 

The  velocity  history  and  drag  versus  Mach  number  data  for  a  typical  fin-stabilized 
(FS)  SFRJ  projectile  are  shown  in  Figures  16  and  17.  The  solid  fuel  in  these  projectiles 
contained  AP.  The  effect  of  ramjet  burn  are  shown  as  nearly  a  constant  Co  between  Mach 
number  3.6  and  4.0  in  Figure  17  and  gradual  retardation  in  velocity  between  0  and  1.6 
seconds  (burn  time)  in  Figure  16.  The  FS-SFRJ  has  a  longer  combustion  chamber  (fuel 
grain)  than  the  spin-stabilized  version.  In  addition,  the  FS-SFRJ  has  nearly  zero  axial 
spin  in  flight. 

The  ramjet  thrust  is  defined  from  this  data  as. 

Thrust  =  (13) 

where  p  is  the  air  density,  V  is  the  velocity  along  the  trajectory,  Srtf  is  the  reference  area 
with  a  reference  diameter  of  75mm,  and  AC'd  is  the  difference  in  drag  coefficient  between 
five  and  inert  SFRJ  projectiles  at  the  same  Mach  number.  Figure  18  shows  the  thrust  as 
determined  from  Equation  13.  The  SFRJ  projectile  generates  about  250  pounds  of  thrust 
during  the  1.6  second  burn  time. 

The  present  steady-state  combustion  model  can  predict  the  SFRJ  burn  time  from  the 
soUd  fuel  regression  rate  and  the  thickness  of  the  solid  fuel  layer, 
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(14) 


burn  time(sec)  = 


ifti(in) 


r{in/sec) 

where  tfa  and  r  ceui  be  average  values  over  the  fuel  grain.  The  drag  coefficient  is  computed 
from  the  change  in  momentum  from  the  injector  (combustion  chamber  inlet)  to  the  nozzle, 


AMomentum  ^  /  /s,  pi  •  hi)dS  -  J  fs,  P2^2(^2  •  n2)dS 

l/2pV^S^  =Cd  =  — ■ - l/2pV^s2 - 

where  subscripts  1  and  2  refer  to  the  injector  and  nozzle,  respectively.  For  the  inert 
(unfueled)  projectile,  for  which  the  internal  flow  is  supersonic  and  an  oblique  shock  wave 
is  attached  to  the  leading  edge  of  the  inlet,  Cd  was  obtained  from  Reference  5.  For  the 
live  projectile,  the  internal  flow  is  subsonic  with  combustion  and  with  an  attached  normal 
shock  at  the  leaxiing  edge  of  the  inlet.  Co  was  obtained  from  Reference  5  for  the  exterior 
and  inlet  and  from  the  present  code  for  the  combustion  chamber. 


The  computed  results  are  shown  with  measurements  in  Figure  18.  The  computed 
bum  time  is  1.9  seconds,  18%  longer  than  the  1.6  seconds  observed.  The  computed  thrust 
is  210  pounds,  16%  smaller  than  the  average  measured  thrust  of  250  pounds.  Although 
the  projectile  decelerates  dvuing  the  burn  time,  the  computation  assumes  a  steady-state 
bum  and  thus  the  launch  velocity  of  4500  ft/sec  is  used.  Since  fuel  oxidizer  (AP)  effects 
are  not  included  in  the  combustion  model  the  computation  overestimates  the  burn  time 
and  imderestimates  the  thrust.  Connected-pipe  tests^^  show  that  the  addition  of  AP  (with 
catalyst)  to  solid  fuel  increases  the  regression  rate  by  about  40%  and  thus  decreases  the 
bxim  time  by  about  40%  in  a  non-spiixning  SFRJ.  Thmst  measurements  are  not  available, 
however  AP  should  increase  the  thmst  generated  by  an  SFRJ  due  to  an  increased  fuel  flow 
rate.  Flight  tests  of  solid  fuels  without  AP  in  SFRJ  projectiles  have  not  been  reported. 
Sources  of  discrepancy  between  measured  and  computed  flight  performance  include:  the 
effects  of  radiation  which  are  excluded  from  the  combustion  model,  vmcertainty  in  the 
thermodynamic  input  quantities,  and  the  validity  of  the  diffusion  flame  model.  Efforts 
should  be  made  to  incorporate  the  effects  of  AP  for  non-spinning  flight  cases. 


VII.  CONCLUSIONS 


A  steady  flow  solution  of  the  Navier-Stokes  equations  is  linked  with  a  diffusion-flame 
combustion  model  for  a  one-step,  irreversible,  mass-controlled  reaction  between  solid  fuel 
Eind  eur.  The  reaction  is  assumed  to  go  to  completion  and  to  produce  products  that  are 
treated  as  a  single  specie.  A  one-equation,  algebraic  turbulence  model  for  combustor  flows 
is  employed.  Solutions  produced  by  this  code  are  applied  to  the  subsonic  flow,  solid  fuel 
combustion  within  a  SFRJ  projectile.  The  results  provide  estimates  of  the  major  features 
of  the  internal  flow  including  flow  recirculation,  flame  position,  exit  plane  flow,  burn  time, 
and  thmst.  The  quality  of  the  results  are  dependent  on  the  accuracy  of  the  thermodynamic 
input  data.  Since  the  code  assumes  a  mass-controlled  reaction,  reaction  rate  data  is  not 
required.  Thus  the  code  can  be  applicable  to  the  wide  range  of  solid  fuels  commonly  used 
in  SFRJ.  Computer  requirements  on  a  Cray  X-MP/48  range  from  .25  to  1.0  CPU  hours  per 
solution,  making  the  code  suitable  for  SFRJ  combustor  design  studies.  This  application 
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of  a  diffusion-flame  combustion  code  serves  as  a  logical  flrst  step  to  detailed  combustion 
modeling  of  reacting  flow  in  SFRJ  projectiles.  The  addition  of  radiation  and  fuel  oxidizer 
(AP)  effects  could  improve  the  utility  of  the  simulation. 

Major  observations  on  the  combustion  in  a  SFRJ  from  computed  results  include: 

1.  The  injector  serves  as  a  flameholder,  in  that  the  flame  is  attached  to  the  injector  and 
extends  downstream  into  the  combustion  cheimber. 

2.  At  any  axial  position  along  the  chamber  the  location  of  the  combustion  flame  is 
indicated  by  a  temperature  maximum,  a  density  minimum,  and  the  mass  fraction  of 
fuel,  oxidizer  and  products  equal  to  0,  0,  and  1  respectively. 

3.  If  the  flame  does  not  reach  the  chamber  axis  at  the  nozzle  throat,  the  flow  at  chamber 
exit  will  consist  of  accelerated  air,  flame,  diluted  and  unbumed  fuel,  and  combustion 
products.  Thus  the  exit  flow  may  not  be  entirely  energized  by  combustion. 
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Figure  1.  SFRJ  Geometry;  Dimensions  in  Millimeters,  Anj  =  43.2  mm,  D^t  =  40.6  mm. 
(I)  Inlet,  (II)  Injector,  (III)  Combustion  Chamber,  (IV)  Fueled  Wail,  (V)  Nozzle. 
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( 1 )  RECIRCULATING  (SEPERATED)  FLOW 

(2)  PRODUCTS  AND  UNCONSUMED  OXIDIZER 

(3)  FLAME,  PRODUCTS  ONLY 

(4)  PRODUCTS  AND  UNCONSUMED  FUEL 


Figure  2.  Schematic  of  Steady-State  Combustion  in  the  Chamber  of  an  SFRJ. 
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Variation  of  Laminar  Reattachment  Length 

with  Reynolds  Number  (Inert  Flow) 
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Figure  3.  Variation  of  Laminar  Reattachment  Length  with  Reynolds  Number,  Inert  Flow, 
Re  =  40  to  200 
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Reattachment  Length  (step  heights) 


Variation  of  Reattachment  Length  with  Spinrate 
in  a  Sudden-Expansion  Combustor  (Inert  Flow,  Re  =  31000) 
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Figure  4.  Variation  of  Reattachment  Length  with  Spinrate,  Inert  Flow.  Re  =  31000. 
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Regression  Rate  (in/sec 


Regression  Rate  Distribution  Along  the  Fueled  Wall 

for  a  SFRJ  Connected-Pipe  Test 


Distance  Fronn  Injector  (in) 


Figure  5.  Regression  Rate  Distribution  Along  the  Fueled  Wall  of  a  SFRJ;  Connected-Pipe 
Test. 
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Figure  6.  Variation  of  Regression  Rate  with  Ma^s  Flow  Rate  in  a  SFRJ;  Connected-Pip^ 
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Figure  8.  Computed  Velocity  Vectors  for  the  Flowfield  Near  the  Nozzle  of  a  SFRJ. 
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Figure  10.  Computed  Streeimwise  Velocity  Profiles. 
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Figure  11.  Computed  Mass  Fraction  Profiles;  X  =  0.136  inches. 
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Figure  12.  Computed  Mciss  Fraction  Profiles;  X  =  2.038  inches. 
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Figure  13.  Computed  Mass  Fraction  Profiles;  X  =  5.942  inches. 


Figure  14.  Computed  Temperature  Profiles. 
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Appendix  A 


In  this  appendix  the  conservation  equations  for  mass,  momentum  and  energy  in  a 
compressible  fluid  flow  are  described.  In  addition,  equations  for  the  conservation  of  chemi¬ 
cal  species  and  various  fluid  properties  axe  derived.  These  equations,  along  with  a  suitable 
model  of  turbulence  viscosity,  form  the  basis  of  the  present  combustion  model.  The  inter¬ 
nal  geometry  of  the  SFRJ  (Figure  1)  promotes  flow  with  substantial  regions  of  low  speed 
flow  (Mach  number  <  1)  with  recirculation  and  swirl  (for  non-zero  axial  spin).  Thus  the 
conservation  equations  axe  elliptic  in  form.  Due  to  the  axisymmetric  nature  of  the  SFRJ 
geometry,  the  equations  axe  written  in  cylindrical  coordinates.  The  velocity  components 
in  this  system  axe  u,u,  and  w  for  the  radial  (r),  azimuthal  (6),  emd  axial  (z)  directions, 
respectively.  Axisymmetric  flow  is  assumed,  thus  all  ^-derivatives  aure  ignored.  However, 
the  azimuthal  velocity  component  (v)  and  the  azimuthal  momentum  equation  are  retained 
so  that  swirling  flows  can  be  treated.  Only  steady-state  flows  axe  considered,  thus  time 
derivatives  (5/ dt)  axe  ignored  as  well. 


For  flows  that  axe  subsonic  and  dominated  by  recirculation  and  rotation,  the  stream 
function-vorticity  form  of  the  governing  equations  has  been  widely  utilized.^®  In  this  form 
an  efiicient  Gauss-Seidel  relaxation  algorithm  is  the  natural  choice  of  solution  procedures. 
The  non-primitive  variables  (stream  function,  ip,  and  vorticity,  u>)  are  defined  as;^® 


w  =  V  X  y  = 
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dz 
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The  conservation  equations  for  mass  and  momentum  are  given  by:^ 
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The  shear  stress  (r )  includes  the  Reynolds  stress  with  an  effective  fluid  viscosity  expressed 
as  the  sum  of  molecular  and  turbulent  viscosities. 


Equation  A4  is  transformed  into  the  non-primitive  variables  using  Equations  Al  and 
A2.  The  resulting  equation  demonstrates  that  the  stream  function,  by  definition,  satisfies 
conservation  of  mass. 


d 

dip 

d 

dip 

dz 

dr 

dr 

dz 

(AS) 


Schlickting,  H.,  Boundary  Layer  Theory,  7th  ed.,  translated  by  J.  Kestin,  McGraw-Hill,  New  York,  1979. 
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Equations  A5  and  A6  are  transformed  into  a  single  vorticity  equation.  First,  the 
definition  of  is  used  in  these  equations  and  yields: 


1  d  /  1  dip  drp 
r  dz  \rp  dz  dr 


dr  \rp  dz  dz 


—  rAin 

r  dz  \rp  dr  dr  j  dr  yrp  dr  dz  j  dz 

The  r-derivative  of  Equation  AlO  is  subtracted  from  the  z-derivative  of  Equation  A9  and 
the  definition  of  vorticity  is  introduced  using  Equation  A3. 
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The  right  hand  side  of  Equation  All  is  expanded  and  the  definition  of  vorticity  (Equation 
A3)  is  again  introduced  to  yield: 
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where, 

VA  =  {dA/dr)f  +  {dA/dz)z 
r  =  unit  vector  in  the  r-direction 

z  =  unit  vector  in  the  z-direction 
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Rearrangement  and  expansion  of  the  gradient  terms  yields  the  final  form  of  the  vorticity 
equation: 


(A14) 


The  azimuthal  momentum  equation  (Equation  A7)  can  be  expressed  in  terms  of  the 
non-primitive  variables  as: 
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Energy  conservation  for  a  compressible  flow  is  expressed  by  the  First  Law  of  Thermo¬ 
dynamics.  The  steady-state  form  of  the  First  Law  states  that  the  net  rate  of  stagnation 
enthalpy  inflow  for  a  control  volume  is  equal  to  the  sum  of  the  shear  work  done  by  the  con¬ 
tents  of  the  control  volume  on  the  surroundings  and  the  heat  transfer  to  the  surroundings.^® 
Turbulence  kinetic  energy  flux  is  also  included  in  the  balance,  while  radiation  flux  is  ne¬ 
glected. 


pVh  +  ^  hjjj  -f  jjfe  -  (ufv  -I-  vfe  -h  wTg) 


=  0 


(A16) 


where, 

h  =  the  stagnation  or  total  enthzdpy  of  the  mixture  [h  =  h  +  /2) 

Jh  =  — (r/i)effCpVT,  the  heat  flux 

Th  =  thermal  conductivity  /  specific  heat  of  the  mixture  =  ^eff/Prandtl  Number 
Cp  =  til®  specific  heat  of  the  mixture  at  constant  pressure 

rrij  =  the  mass  fraction  of  species  j 
T  —  temperature 

hj  =  the  heat  of  formation  of  species  j 

Jj  =  — (rj)effVmj,  the  mass  flux  of  species  j  per  unit  area  per  unit  time 
Fj  =  mEiss  diffusion  coefficient  x  mixture  density  =  ^eff /Schmidt  Number 
Jk  =  —  (rfc)efrVA:,  the  flux  of  turbulence  energy 

Tk  =  turbulence  kinetic  energy  diffusion  coeff.  x  mixture  density  =  /^eff/Schmidt  Number 
k  =  kinetic  energy  of  the  turbulent  fluid  motion  {V'^/2  is  the  KE  of  the  mean  motion) 


For  turbulent  flows,  the  subscript  “eff”,  which  stands  for  “effective”,  serves  as  a  reminder 
that  the  turbulent  fluid  properties  axe  being  treated  as  kinematic  properties  with  a  set  of 
augmented  transport  coefficients.  The  diffusion  coefficients,  F;,,  Fj,  and  Ft  are  defined  in 
terms  of  either  the  Prandtl  Number  (cpp/n)  or  the  Schmidt  Number  {fi/{pDj)).  These 
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nondimensional  numbers  represent  the  ratio  of  momentum  transport  to  thermal  conduc¬ 
tivity  (k)  and  the  ratio  of  viscous  momentum  transport  to  diffusive  m2iss  transport  {Dj  = 
diffusion  coeflBcient),  respectively. 

The  energy  equation,  Equation  A16,  is  transformed  into  stre8im  function-vorticity 
form.  First,  the  terms  involving  shear  stress  in  Equation  A16  can  be  expanded  to  yield: 
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If  the  terms  in  square  brackets  of  Equations  A17-A19  are  small  and  can  be  set  to  zero 
(which  is  exactly  the  case  on  a  solid  wall  boimdary),  the  shear  stress  term  in  Equation 
A16  can  be  more  easily  written  as: 


(uTr  +  VTg  -I-  WT^)  =  /XeffV(V^/2)  =  /UeffV 


-f 


(A20) 


This  assumption  will  significantly  ease  subsequent  manipulations,  and  indeed  Equation 
A20  serves  as  a  good  approximation  for  most  practical  applications. “-25  The  equation  for 
energy  conservation  is  then  given  by: 
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Equation  A21  can  be  expressed  in  the  non-primitive  variables  as: 
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where. 


Oj  =  A^eff/Fj,  Schmidt  Number 
<^h  =  Peti/Thi  Prandtl  Number 

<^k  =  A^eff/Fjt,  Schmidt  Number 


For  a  steady-state  flow  the  conservation  of  chemical  species  states  that  the  sum  of 
the  net  inflow  of  species  j  by  convection,  the  net  inflow  by  diffusion,  and  the  creation  of 
species  j  by  chemical  reaction  (iZj,  per  imit  volume)  within  the  control  volume  is  equal  to 
zero. 
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Or  in  non-primitive  variable  form: 
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The  governing  equations  for  stream-fimction  (Equ.  A3),  vorticity  (Equ.  A14),  az¬ 
imuthal  velocity  (Equ.  Al5),  stagnation  enthalpy  (Equ.  A22)  and  species  mass  (Equ. 
A24)  can  be  expressed  in  a  general  form  as  follows: 


d  (  ,  dtl)\  d  (  ,  dtl>\ 

d 

1  ^ 

d 

1  ^ 

Tz  vd^J 

dz 

M) 

dr 

b^r—  {c^4>) 

-h  rd^  =  0 


(A25) 


with. 


II 

-o- 

=  0, 

Cp  =  1, 

1 

II 

II 

-o- 

=  r'^, 

It 

Cp  =  fCfg, 

dp  — 

II 

dtj)  1, 

bp  = 

Cp  =  l/r^, 

=  0 

II 

-o- 

=  1, 

bp  =  (r/,)eff  , 

C0  =  1, 

dp  =  Qh  —  Sh 

<f>  =  rrij, 

Clip  =  1, 

bp  =  (rj)ejr. 

Cp  =  1, 

dp  —  JRj 

where, 


The  pressure  can  be  recovered  from  the  non-primitive  variables  after  a  converged 
solution  of  Equation  A25  is  achieved  or  after  each  iteration,  if  pressure  variations  are 
expected  to  have  a  significant  effect  on  density.  The  radial  and  axial  momentum  equations 
(Equ.  A5  and  A6)  can  be  rearranged  to  yield: 

dp  d 

—  =  Pi(r,z)  and  -£  =  P2{r,z)  (A26) 
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Then  along  any  path  from  point  A  to  point  B  in  the  flowfield,  the  pressure  is  given  by; 

ffl 


Pb-Pa=  f  (Pidz  +  P2dr) 
J  A 


(A27) 


Since  p  is  a  scalar,  pB  —  Pa  should  be  path  independent  and  can  serve  as  a  consistency 
check  on  the  converged  solution.  In  most  cases  a  pressure  difference  is  desired  to  form  the 
pressure  coeflBcient  at  a  point.  However,  if  the  pressure  at  point  B  is  required,  the  pressure 
at  point  A  is  assigned  to  a  known  inlet  pressure  value,  and  integration  proceeds  from  the 
inlet  to  point  B. 


A  thermodynamic  equation  that  relates  the  mixture  density,  p,  and  the  mixture  tem¬ 
perature,  T,  to  the  stagnation  enthalpy,  is  given  by. 


dh  =  ^  (cpjirij)  dT  -t-  hjdrrij  +  ^  +  dk 


(A28) 


where,  d{V^)/2  and  dk  represent  the  kinetic  energy  from  the  mean  flow  and  the  fluctuating 
tmrbulent  flow,  respectively.  The  heats  of  formation,  hj,  must  be  specified  from  a  knowledge 
of  the  fuel  chemistry.  Assuming  that  hj  is  evaluated  at  To  =  0  and  that  the  specific  heats 
are  independent  of  T  (i.e.  calorically  perfect  gases),  then  Equation  A28  can  be  integrated 
to  yield: 

h  =  T  ^  Cpjrrij  -I-  J]  hjinj  ■¥  —  +  k  (A29) 

3  i  ^ 

The  csdorically  perfect  gas  assumption  is  made  to  ease  the  computational  effort  and  is  not 
unusual,^^  but  it  could  be  removed  if  the  temperature  dependence  of  Cp  for  the  reactants 
and  products  were  known.  The  temperature  of  the  mixture  is  then  expressed  as; 

j,  _  -k-J2i  hjTrij 

Ej  Cpjrrij 


The  specific  volume  of  an  individual 
defined  as; 


species  and  the  specific 


PjMj 


volume  of  the  mixture  are 

(A31) 


15  =  I3jmj 


diTruj 

PjMj 


(A32) 


where  3^  is  the  universal  gas  consteint  and  M,  is  the  molecular  weight  of  species  j  (note, 
^  Hj  ^jTTij  since  all  j  components  occupy  the  seime  volume).  The  mixture  density  is 
defined  using  p  =  H,jPj  =  ^'PHj{^j/Mj)/l5  and  Equation  A30. 


P  y 


1  p/3?r 


p/3i 


Uj 


'■-f 


m. 


Ej  Cpjirij 

E171 1 

3  M, 


(A33) 
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Appendix  B 


Finite  difference  formulae  axe  employed  to  transform  Equation  A25  into  an  algebraic 
equation  in  <f>  for  a  single  computational  cell.  The  restJt  for  all  cells  in  the  computational 
mesh  is  a  system  of  equations  that  can  be  solved  eflBciently  using  a  Gauss-Seidel  relaxation 
scheme.  As  outlined  below,  the  computational  stability  of  the  scheme  is  insured  since 
upwind  differencing  is  employed  for  the  convection  terms.  Under-relaxation  is  used  to 
accelerate  convergence  aind  enhance  stability  when  the  magnitude  of  the  source  terms 
(especially  with  combustion)  is  large. 


Figure  Bl  shows  the  computational  mesh  in  the  vicinity  of  a  typical  mesh  point,  P. 
The  integration  of  Equation  A25  is  carried  out  on  a  rectangle  with  comer  points  denoted 
ne,se,sw,nw. 


(Bl) 


Equation  Bl  can  be  integrated  if  is  constant  which  is  the  case  for  all  (j)  except  when 
<l>  =  u/r  (i.e.  =  r^).  Assuming  that  is  a  suitable  average  value  for  over  the 

rectangle  around  point  P,  then  Equation  Bl  is  rewritten  as: 


where  the  subscripts  n,e,s,w  when  applied  to  a  spatial  derivative  indicate  that  the  derivative 
is  evaluated  along  the  corresponding  side  of  the  rectangle.  The  terms  in  Equation  B2  are 
due  to  convection,  diffusion  in  the  r-direction,  diffusion  in  the  z-direction  and  the  source 
term,  respectively.  Each  of  these  terms  is  treated  in  a  separate  way. 


Consider  the  first  integral  for  the  convection  terms  in  Equation  B2: 


The  average  value  for  <f>e  is  defined  by: 


(U)^dr  [a^p  (V»ne  -  t^se)] 


(B3) 


(B4) 
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such  that  the  integral  Ic  is  written  as: 

Ic  =  ii’ne  -  V’se)  (B5) 

The  value  of  is  chosen  as  the  value  of  <j>  “upwind”  of  the  e-face  of  the  rectangle.  This 
is  the  upwind  differencing  method  since  the  upwind  direction  depends  upon  the  direction 
of  the  local  velocity  vector.  Thus  if  (V'ne  ~  ^se)  —  {drl?/dr\  >  0  then  <f>c  =  and  if 
(V’ne  -  V’se)  —  (dtp/dr)^  <  0  then  =  ^E-  Thus: 

Ic  ~  a^p^^E  [(V’ne  -  V’se)  “  |  V’ne  “  V'xs  1]  /2  +  a^p<f>P  [(^nt  -  V’se)  +  |  V’ne  “  V’se  I]  /2  (B6) 

where  the  value  of  rp  at  the  rectangle  corner  points  is  a  average  of  the  values  of  t/»  at  the 
neighboring  mesh  points  (i.e.  V’se  —  [V’SE  +  V’e  +  V’P  +  V’s]/4).  Treating  the  other  integrals 
of  the  convection  terms  in  the  same  fashion  we  obtain: 

-^convection  =  (<^p  —  <^e)  +  Aw  (</>p  —  <i>w)  +  An  {<I>p  —  <i>N)  +  As  {<f>p  —  4>s)  (B7) 


where, 

Ae  =  <J^p[(V’*e  V^ne)")"  |  V’se  V’ne  |]/2 

Aw  =  ®0p  [(V’nw  V’sw)"!"  I  V’nw  0sw  |]/2 

An  =  ^0P  [(V’ne  V’nw)"!"  |  V’ne  V’nw  11/2 

As  =  ^<pp[{'^aw  V’se)"!"  I  V’se  |]/2 


Note  that  the  coefficients  Ae,  Aw,  An,  and  As  axe  positive  or  zero. 
Consider  the  first  integral  for  the  diffusion  terms  in  Equation  B2: 

d 


Id  =  J  [^6^1 


dz 


M) 


dr 


(BS) 


Assuming  that  (pe  -|-  rp)/2,  and  that  6^,  ~  (6^^  -t-  b^p)/2  and  using  a  central  finite 
difference  formula  for  {d{c^<t>)/dz)c  ~  -  c<t>p<f>p)f{zE  -  zp),  then; 


A  b,f)p 

>E  +  rp] 

-  C<t>p 

>N  —  rs' 

2 

2  J 

ze 

-  Zy> 

2 

(B9) 


Treating  the  other  integrals  of  the  diffusion  terms  in  the  same  fashion  we  obtain: 

■^diffusion  =  -Be  [c^e'^e  —  c^p  </»?]  + Bw  [c0w‘^w  —  c^p<^p] -H  Bn  —  c^p'ji’p] 

+  Bs  [c,^s<^s  -  c^p'Ap]  (BIO) 


where, 

Be  =  [(fe^E  +  ^'^p)/8][(^n  -  rs)/(^E  -  2p)][rE  +  rp] 
Bw  s  [{b^y,  +  fe0p)/8][(rN  -  rs)/(^P  -  ^w)][rw  +  rp] 
Bn  =  +  ^'^p)/8][(^e  -  ^w)/(rN  -  rp)][rN  rp] 

Bs  =  [(6^s  +  ^■^p)/8][(^e  -  zw)/{rp  -  rs)][rs  -t-  rp] 
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Note  that  the  coefficients  Be,  Bw,  Bn,  and  Bs  are  positive. 
Consider  the  integral  for  the  source  term  in  Equation  B2: 

rrn  fZ* 

/source  =  /  /  rd^dzdr  =  0 

Jzyt 


(BU) 


Assuming  that  d^  is  uniform  over  the  area  of  integration  (rectangle)  and  equal  to  dp  and 
that  rp  is  a  close  approximation  to  the  space-averaged  value  of  r  over  the  area  of  integration 
then; 


d^pTp 


1 

tn  —  rs' 

L  2  J 

[  2 

(B12) 


When  the  expression  for  d^  involves  spatial  derivatives  (e.g.  when  (f>  =  uj/r  and  4>  =  ft), 
central  finite  difference  expressions  are  used.  For  example,  d^  =  {d(j>/dr)p  is  expressed  on 
a  non-uniform  grid  as: 


dr  , 


(</>N  -  <?^p)(rp  -  rs)/(rN  -  rp)  +  {(j>p  -  «?^s)(rN  -  rp)/(rp  -  rs) 


(B13) 


'  p  tn  —  rs 

which  for  a  uniform  grid  reduces  to  (<^n  —  <^s)/(rN  —  rs). 

The  complete  finite  difference  form  of  Equation  B2  for  a  single  grid  point  is  written 


as: 


Ae  (<f>p  —  <^e)  +  Aw  {<f>p  —  </>w)  +  An  (<^p  —  4>n)  +  As  (<^p  —  (^s)  — 

Be  -  c^4>p)  ~  Bw  {c4,^<i>w  —  — 

Bti  -  c*p<^p)  -  Bs  (c,*s^s  -  Cflip'i^p)  +  d^pllp  =  0  (B14) 

which  is  expressed  as  a  successive-substitution  formula  by  moving  <j>p  to  the  left-hand-side 
of  Equation  B14: 

<j>P  =  Ce<^e  +  Cw<f>w  +  C'n<^n  +  Cs(f>s  +  D  (Bl5) 


where, 

Ce  =  (Ae  -I-  Bec^pj)/ X;ab 

Cw  =  (Aw  +  Bwc^w)/  Eab 

Cn  =  (An  -I-  Bnc,^j,  )/  XIab 

Cs  =  (As  -I-  BsC0s)/  J^ab 

B  =  ~<^(*pnp/x3AB 

XIab  =  Ae  +  Aw  +  An  -I-  As  -I-  c^^(Be  +  Bw  +  Bn  -I-  Bs) 

Ae  =  a^p  [{ipsE  +  V’s  -  0ne  -  '0n)+  |  ^se  +  V’s  -  t^NE  -  0n  |]/8 

Aw  =  a^p  [(t^NW  +  0N  —  V’SW  -  V’s)+  1  V’NW  +  0N  —  V’sw  -  '0s  |]/8 

An  =  [(0NE  +  0E  -  0NW  -  0w)+  |  0ne  +  0e  -  0nw  -  0w  |]/8 

As  =  [(0SW  +  0w  -  0SE  -  0e)+  I  0SW  +  0w  -  0SE  -  0E  |]/8 

lip  =  rp(zE  —  2w)(rN  —  rs)/4 

Be  =  [(ft^E  +  &^p)/8][(rN  -  rs)/(2E  -  2p)][rE  +  rp] 

Bw  =  [(ft,^vv  +  ^<<>p)/8][(rN  -  rs)/(2p  -  ^w)][rw  +  rp] 

Bn  =  [{b^  +  ^^p)/8][(2e  —  ■^w)/(rN  —  rp)][rN  +  rp] 

Bs  =  [(ft.^s  +  hv)l^][{zE  -  zw)/{rp  -  rs)][rs  +  rp] 
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Equation  B15  represents  a  successive-substitution  formula  for  <f>p  at  each  grid  point 
on  the  computational  mesh.  The  resulting  system  of  equations  for  the  entire  mesh  is  solved 
using  a  Gauss-Seidel  relaxation  scheme.  Each  iteration  cycle  is  made  up  of  N  sub-cycles, 
where  N  is  the  number  of  equations  (see  Equation  A25)  being  considered  (N  must  be  at 
least  2  since  the  equations  for  <j)  =  uj/r  and  (j>  ■=  x]}  are  the  minimum  required  to  define 
the  flow).  In  each  sub-cycle,  each  grid  point  is  scanned  row  by  row  and  a  single  variable 
is  updated  using  Equation  B15.  The  variables  u/r  and  t/>  axe  updated  in  order,  followed 
by  ^1  other  variables.  When  all  sub-cycles  axe  completed,  a  new  iteration  cycle  is  started 
in  which  the  values  of  the  variables  from  the  latest  iteration  axe  immediately  used.  This 
is  consistent  with  the  Gauss-Seidel  methodology.^®  This  procedure  is  repeated  until 
convergence  is  achieved.  Convergence  is  usually  indicated  by  satisfaction  of  one  or  both  of 
the  following  criteria:* 

<  Ai  (B16) 

max 

<  Aj  (B17) 

max 

where  Ai  and  Aj  axe  small  (Ai  <  .005  and  A2  <  Ai/lOO)  and  superscript  k  represents  values 
at  the  latest  iteration  cycle. 

For  the  successive-substitution  formula  (Equation  Bl5)  to  be  stable  and  convergent, 
three  criteria  must  be  met.^® 


1.  The  sum  of  the  moduli  of  Ce,  Gw,  Gn,  and  Gs  must  be  less  than  or  equal  to  1  at 
every  grid  point. 

2.  The  sum  of  the  moduli  of  Ge,  Gw,  Gn,  and  Gs  must  be  less  than  1  for  at  least  one 
grid  point. 

3.  The  coefficients  (G  and  D)  of  <f)  in  Equation  B15  must  not  vary  greatly  from  one  cycle 
to  the  next. 

The  coefficients,  G,  axe  comprised  of  coeflicients,  A  and  B,  which  are  greater  than  zero 
(a  direct  consequence  of  upwind  differencing  used  for  the  convection  terms,  see  Appendix 
B).  This  implies  that  the  G  coefficients  are  greater  than  zero  and  that  ^  G  =  1  provided 
the  C0  coefficients  axe  equal.  In  fact  \  for  all  definitions  of  (t>  except  when  <i>  ■=  uj/r 
and  4>  =  TV  (Equation  A25).  However,  the  relaxation  scheme  may  be  stable  and  converge 
if  the  variation  in  between  neighboring  grid  points  is  small  and  thus  G  ~  1.  Item  1 
of  the  stability  criteria  will,  therefore,  be  satisfied. 

Along  the  boundaxy  of  the  computational  domain,  there  will  be  at  least  one  case  where 
(j>P  =  constant.  Thus  the  G  coefficients  of  Equation  B15  will  be  zero  and  so  ^  G  <  1  for  the 
grid  points  along  the  boundary.  Item  2  of  the  stability  criteria  will,  therefore,  be  satisfied. 

Satisfax:tion  of  the  third  criteria  cannot  be  shown  formally.  The  computational  grid 
must  be  appropriately  defined  with  clustering  in  regions  of  large  (j)  gradients.  The  initial 
and  boundaxy  conditions  must  be  consistent  with  the  mathematical  nature  of  the  governing 

**  CamaAan,  B.,  L%ihtr,  H.A.,  and  WilkcM,  J.O.,  Apylied  Numerical  Mtihodt,  John  Wiley  and  Sont,  New  York,  1969. 
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partial  differential  equations.  In  addition,  the  variation  in  the  source  terms  must  not 
be  stronger  than  the  damping  of  the  viscous  terms  in  the  equations.  If  divergence  is 
found  to  occur,  then  under-relaxation  techniques  can  be  employed  to  stabilize  the  iteration 
scheme.  In  contrast,  over- relaxation  techniques  tend  to  accelerate  the  convergence  of  a 
stable  scheme,  but  they  are  not  recommended  for  equations  with  strong  source  terms. 
Under-relaxation  is  defined  for  iteration  number  k  by  +  (1  ~  where 

0  <  Our  <  1  and  0^  =  1  for  no  under-relaxation.  For  the  current  governing  equations  and 
combustion  model,  under-relaxation  was  not  required  for  convergence  to  a  stable  solution. 
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Figure  Bl.  Portion  of  the  Computational  Grid.  Solid  circles  are  grid  points.  Open  circles 
axe  cell  centroids.  Dotted  lines  enclose  the  rectanguleir  area  of  integration. 
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List  of  Symbols 


Our  under-relaocation  factor 

Cp  specific  heat  at  constant  pressure 

Cd  drag  coefiicient  (Equation  15) 

D  maximum  diameter  of  the  combustion  chamber 
/  mixture  fraction  (Equations  5-7) 

g  acceleration  due  to  gravity 

h  heat  of  formation 

hy  heat  of  vaporization 

h  total  enthalpy  (Equation  A29) 

J  flow  rate  or  flux 

k  turbulence  kinetic  energy 

L  length  of  the  combustion  chamber 

m  mass  fraction 

M  molecular  weight 

n  stoichiometric  air/fuel  mass  ratio 

n  unit  vector  in  the  normal  direction 
p  static  pressure 

Pr  Prandtl  Number  (see  also  <7^) 

q  heat  flux 

3?  universal  gas  constant 

r  radial  direction 

r  solid  fuel  regression  rate  (Equation  12) 

R  reaction  rate  per  unit  volume 

Re  Reynolds  Number 

S  area 

Sc  Schmidt  Number  (see  also  cTj  and  cr^) 
t  time 

tfu  thickness  of  the  solid  fuel  layer 

T  temperature  (Equation  A30) 

u  velocity  component  in  the  axial  direction 

V  velocity  component  in  the  azimuthal  direction  (swirl  velocity) 

V  magnitude  of  the  local  velocity  vector 

A 

V  ur  +  vO  +  wz 

w  velocity  component  in  the  radial  direction 
z  axial  direction 

Greek  Symbols 

r  diffusion  coefficient 

e  turbulence  dissipation  rate 

9  azimuthal  direction 

K  thermal  conductivity 

fi  molecular  viscosity  (Equation  10-11) 
p  density  (Equation  A33) 
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ch 

Pca/Thi  Prandtl  Number 

/^eff/Pj,  Schmidt  Number 

Ok 

Peal^ki  Schmidt  Number 

T 

shear  stress 

<i> 

general  flow  variable  (Equation  A25) 

stream  function  (Equation  A1-A2) 

(jj 

vorticity  (Equation  A3) 

15 

speciflc  volume  (Equation  A31-A32) 

Superscripts 

* 

unit  vector 

total  or  stagnation 

rate 

Subscripts 

eff 

effective 

fo 

fuel-oxidizer 

fuel-products 

fu 

fuel 

h 

enthalpy 

j 

j-th  mixture  component 

k 

turbulence  kinetic  energy 

op 

oxidizer-products 

ox 

oxidizer 

p 

constant  pressure 

pr 

products 

r 

radial  component 

ref 

reference 

St 

stoichiometric 

t 

turbulent 

vap 

vaporization 

WElll 

on  the  wedl 

z 

axial  component 

9 

azimuthal  component 
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